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Systems whose underlying classical dynamics are chaotic exhibit signatures of the chaos in their 
quantum mechanics. We investigate the possibility of using time-dependent density functional theory 
(TDDFT) to study the case when chaos is induced by electron-interaction alone. Nearest-neighbour 
level-spacing statistics are in principle exactly and directly accessible from TDDFT. We discuss 
how the TDDFT linear response procedure can reveal the mechanism of chaos induced by electron- 
interaction alone. A simple model of a two-electron quantum dot highlights the necessity to go 
beyond the adiabatic approximation in TDDFT. 
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I. INTRODUCTION 

The study of the quantum mechanics of systems whose 
underlying classical dynamics are chaotic, has revealed 
many intriguing features. It is now well-established that 
underlying classical chaos dramatically manifests itself in 
certain quantum signatures: spectral fluctuations, local- 
ization properties of wavefunctions, inverse participation 
ratios, and extreme sensitivity to tiny variations in con- 
trol parameters p}. Perhaps the most widely studied of 
these is the spectral statistics: In particular, the nearest- 
neighbour spacing distributions (NNS) for intcgrable sys- 
tems generically display Poissonian statistics (level clus- 
tering), while displaying Wigner-Dyson statistics (level 
repulsion) for chaotic systems. Such correspondence be- 
tween classical chaos and NNS statistics was conjectured 
by Bohigas, Giannoni, and Schmit Q, and while coun- 
terexamples exist [H, Hj , the association is so general that 
Wigner-Dyson statistics are often viewed as signatures 
of underlying classical chaos. This property has been 
used to study the transition from integrability to chaos 
as a parameter in the confining potential of the system 
is varied |5j. Many of these studies consider one parti- 
cle, or, in quantum dots, non-interacting electrons; the 
chaos arises in such cases due to the shape of the dot 
confining potential. What is less well-understood is what 
happens when electron-interaction is turned on. Suppose 
that interacting electrons are placed in a potential where 
the single-particle NNS is Poissonian; is the chaos in- 
duced by the Coulomb interaction enough to transform 
the NNS statistics to a Wigner-Dyson distribution? To 
what extent does the level repulsion kick in? Is the pic- 
ture qualitatively different with few electrons compared 
to many electrons? 

That electron-interaction alone can induce chaos is cer- 
tainly evident from the very early days of quantum me- 
chanics, impeding Bohr in 1913 from successfully quan- 
tizing the Helium atom [6| (task only completed semiclas- 



sically in 1991 0); for a beautiful review of the theory of 
two-electron atoms, see ref.Q). Putting aside concerns 
regarding the validity of random-matrix theory when 2- 
body interactions are present we note that the first 
application of random-matrix theory outside of nuclear 
physics was to a series of complex atoms (ill fl2j , where 
the spin-orbit interaction provided the crucial ingredi- 
ent in yielding Wigner-Dyson statistics. Statistics of a 
different series of complex atoms based on experimen- 
tal data [HI], as well as theoretical models [l4[ support 
the finding that highly excited states of complex atoms 
tend to display Wigner-Dyson statistics. There have been 
several studies in molecules (eg. (lo| ) where the coupling 
of electronic excitations with nuclear vibrational and ro- 
tational excitations provide the complexity. Generally 
experimental data, with many levels of the same symme- 
try is needed. The chaos in these examples is understood 
to arise not from Coulomb electron-interaction alone but 
rather from its coupling to other degrees of freedom. 

For simple atoms, there have been very few calcu- 
lations of the level-spacing statistics, perhaps because 
of the challenges involved in gathering enough levels of 
doubly-excited resonances, either from theory or exper- 
irnent, which appear to be crucial for this effect (see 
Ref. [Ty] for a calculation of level-repulsion in helium). 
However deviations from Poissonian statistics have been 
clearly identified in non-hydrogenic Rydberg atoms in 
either magnetic or electric fields [13, E, 1 19L [20| . The 
hydrogen atom at the corresponding parameter regime 
is integrable, but scattering off the ionic core in a non- 
hydrogenic Rydberg atom creates fundamentally differ- 
ent dynamics, with chaotic trajectories depending on the 
value of the quantum defect. In models of this effect, 
NNS have been shown to display level-repulsion, follow- 
ing a distribution intermediate between Poissonian and 
Wigner-Dyson [l7l |. 

There are also many-electron solid-state examples, 
where electron-interaction has shown to lead to the 
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transition between Poissonian and Wigner-Dyson statis- 
tics pll . [22l . I23I [2^ |. and several of these works iden- 
tify the driving parameter for this transition. Difficulties 
with interpreting experimental data make the idea of a 
theoretical calculation of the level statistics attractive. 
At the same time, because the solution of the interact- 
ing many-electron problem grows exponentially with the 
number of electrons, typically a model is used for the in- 
teraction, eg. a two-body random interaction model was 



used in Ref. [2l|. A recent calculation of just two elec- 
trons in a quantum dot [25| , explicitly demonstrates that 
Coulomb electron-interaction alone can transform the 
Poissonian statistics of the non-interacting system into 
Wigner-Dyson. Generally in few-electron quantum dot 
studies of transport (26j |. random constant-interaction 
models [13] are often used, however Ref. [25| suggests 
this may not be a good description of the real interacting 
electronic system, as constant-interaction models retain 
the Poissonian statistics of the integrable non- interacting 
system. 

In order to better understand the mechanisms of 
interaction- induced chaos, it would be desirable to use 
a method that captures electron correlation reliably and 
efficiently (given the thousands of excited states needed 
in the calculation), scales well with the number of elec- 
trons, and from whose procedure one can glean aspects 
of the mechanism that brings about the chaos. 

Time-dependent density functional theory (TDDFT) 
is the leading candidate for such a method. This has 
become the method of choice for the calculation of ex- 
citations and response properties of interacting elec- 
tronic systems, because of its scalability: typically the 
accuracy is comparable to sophisticated wavefunction 
methods such as Complete Active Space Self- Consistent 
Field (CASSCF), while implementations follow a far 
cheaper scaling with system size, comparable to time- 
dependent Hartree-Fock (TDHF). The theory in prin- 
ciple yields exact excitations, but approximations for 
exchange-correlation effects are needed in practise. Typ- 
ically, excitation energies are given to within a few tenths 
of an eV, although there are notorious exceptions [Iol.[28|. 
There has been a tremendous drive in recent years to de- 
velop and improve the currently-available functionals. 

In this paper we will explore the possibilities of us- 
ing TDDFT for investigating quantum chaos induced 
by electron-interaction alone. This is a new area for 
TDDFT. The idea is to use TDDFT to study the tran- 
sition from clustering to repulsion statistics in a given 
system as the electron interaction is turned on. For 
the case of chaotic quantum dots, ground-state DFT has 
been used within the local spin density approximation 
to study the statistics of ground-state spin and spac- 
ing between conductance peaks in the Coulomb block- 
ade regime (addtion spectra) [2!| [3(| , but to our knowl- 
edge, the excited-state statistics of isolated dots with 
fixed number of electrons have not been studied via den- 
sity functional methods. We will show in Section HT1 how, 
by monitoring the evolution of the spectral statistics at 



different stages of the calculation, the TDDFT linear re- 
sponse framework can also shed light on the mechanism 
of interaction-induced chaos. 

A most essential question is whether the present-day 
functionals are good enough to perform these tasks. In 
Sec. IIIII we will give an example of a case where they 
are not. This example highlights an important challenge 
that must be overcome for TDDFT to be used in these 
studies. 

Aside from fundamental interest, it is important to 
characterize chaotic versus integrable dynamics for appli- 
cations of technological interest, such as transport across 
quantum dots, highly excited atoms and molecules in ex- 
ternal fields, quantum control and manipulation in ex- 
ternal fields, or the engineering of quantum computer 
hardware that maximizes the fidelity of quantum com- 
putations (3lj . 



II. LEVEL STATISTICS FROM 
TIME-DEPENDENT DENSITY FUNCTIONAL 
THEORY 

TDDFT [32] has become the method of choice to calcu- 
late a variety of response properties of molecules, clusters, 
and solids, in the presence or absence of external time- 
dependent fields [28| . The Runge-Gross theorem provides 
a rigorously exact foundation for the theory: this states 
that given the initial state of the interacting electronic 
system, all observables of the system can be extracted in 
principle from just the time-evolving density. 

Most applications in chemistry and solid-state physics 
currently fall in the linear response regime, where 
TDDFT yields predictions for the optical spectra, i.e. 
the frequency and intensity of electronic excitations in 
response to electric fields [H, H3] . 

We now explain in some detail the standard compu- 
tational procedure. Consider a time-independent N- 
electron Hamiltonian (atomic units will be used through- 
out) 



H = T + V cl 



drn(r> ext (r) 



(1) 



whose energy spectrum {Ek} we want to calculate. In 
Eq.Q f = -|X)iVf stands for the TV-electron ki- 
netic energy operator, y oc = | YIum l r » — r jl _1 f° r the 
electron-electron repulsion, u ex t(r) for the external po- 
tential due to the nuclei, or applied static fields, and 
n(r) = J2iLi ^( r — *i) i s the density operator. The first 
step of a linear-response TDDFT calculation involves the 
self-consistent solution of the ground-state Kohn-Sham 
(KS) equations 



1. 



v s [n]{r) 4>i(r) = e»0,-(r) 



(2) 



Here v s [n] (r) is the KS potential defined such that N non- 
interacting electrons in v s [n] (r) have the same ground- 
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state density n(r) = (>I'o|rj( r )| , I'o) as the original inter- 
acting system of ground state l^o)- hi the KS scheme, 
the density is obtained from the N occupied KS orbitals 
as n(r) = ^ i occ \<fii(r)\ 2 . The KS potential is written as 
the sum of three pieces: 

u s [n] (r) = u ext (r) + v H [n] (r) + v xc [n] (r) , (3) 

where u H M( r ) = / d 3 r'n(r')/\r — r'\ is the Hartree po- 
tential, and u X cM( r ) the exchange-correlation poten- 
tial. This is the functional derivative (with respect to 
the density) of the exchange-correlation energy func- 
tional E xc [n), evaluated at the ground-state density: 
ti xc [n](r) = SE XG [n]/dn(r)\ . The functional E xc [n] is 
the only quantity that needs to be approximated in or- 
der to get the ground-state density and energy; it is for- 
tunately amenable to local approximations [33 ] in many 
situations. Knowledge of u xc (r) is sufficient to calculate 
the ground-state energy Eq, but the excited-state ener- 
gies {Ek}, k > 0, are not accessible from this first step. 
The occupied orbital energies and KS orbitals 4>i(r), 
along with the unoccupied orbital energies and orbitals, 
are used as basic ingredients for the second step of the 
calculation. The aim of this second step is to correct the 
unphysical KS excitations towards the correct ones. Here 
one solves for eigenvalues and eigenvectors of the matrix: 

^99' = <W + 4 v / £y^7 [q | /hxc M I q'} . (4) 

The square-root of the eigenvalues correspond to the ex- 
citation frequencies E^ — Eq. The double index q = (i, a) 
represents a single excitation from occupied KS orbital 
4>i{r) to unoccupied KS orbital 4> a {r), w q is the difference 
between occupied and unoccupied KS orbital energies, 
u q — e a - ei, and 

[gl/nxoMk'] = / drdr'<f>*(T)<l> a (r) 

x/ HXC (r,r';c)^(r')^(r') . (5) 

The Hartree-exchange-correlation kernel, 

/HxoH( r , r '; w), is the central quantity of linear-response 
TDDFT. In the time-domain, f HXC [n](rt,r't') is the 
functional derivative of the time- dependent Hartree-plus- 
exchange-correlation potential v H [n t ](r,t) + v xc [n t ](r,t) 
with respect to the time-dependent density n t (r't'), 
evaluated at the ground-state density n(r')(see Sec. 4. 3. 2 
of [36[ for more details). Eqs. (g]) and (0 are obtained 
from the linear response limit of the full TDDFT equa- 
tions, in which N electrons evolve in the time-dependent 
KS potential 

u s [n t ] (r, t) = v cxt (r, t) + v H [n t ] (r , t) + v xc [n t ] (r,t) . (6) 

Here, nt = n(r, t) is the time-dependent density of the N 
KS electrons, which is the same as that of the interacting 
system [32j]. 

To summarize, there arc two stages: (1 st ) The 
ground-state Kohn-Sham (KS) equations are solved 
to self-consistency; this requires an approximation to 



the ground-state exchange-correlation energy functional. 
Even if the (unknown) "exact" functional were used here, 
the excitations of the ground-state KS potential (i.e dif- 
ferences between occupied and unoccupied KS orbital en- 
ergies) could only be regarded as zeroth-order approxi- 
mations to the true excitations of the system. (2 nd ) The 
KS frequencies are corrected via / HX c to become the true 
excitations of the many-body system [33| . 

With rare exceptions, an adiabatic approximation 
(ATDDFT) is employed in this second step: the 
exchange-correlation potential at time t is approximated 
by that of a ground-state of the instantaneous density 
n(rt). This means that the Hartree-exchange-correlation 
kernel has no frequency-dependence in ATDDFT. An im- 
portant consequence for the purposes of this paper is that 
states of double (or higher multiple) excitations cannot 
be captured within ATDDFT [12, HH ; such excitations 
require a strongly frequency-dependent kernel. In the 
simplest (and common) adiabatic local-density approxi- 
mation (ALDA), the approximation is local in space as 
well as time. 

Notwithstanding the resounding success of TDDFT 
within existing approximations, different instances where 
the approximations fail have been identified. Charge- 
transfer excitations at long-range, conical intersections, 
states of multiple-excitation character, polarizabilities of 
long-chain polymers, Rydberg excitations, lie among the 
challenges in the linear response regime; non-scqucntial 
photoionization and quantum control applications are 
important challenges for approximations in the strong- 
field regime. Given the success of TDDFT for the vast 
majority of problems, and the fact that it scales in 
a reasonable way with system size while incorporating 
electron-correlation effects, there is a tremendous drive in 
the recent literature to understand and improve TDDFT 
approximations. 

Elucidating the mechanism of interaction- 
induced chaos: type-KS vs. type-/ HXC chaos 

Consider now a system of N interacting electrons with 
Wigner-Dyson NNS statistics, and focus on the case 
where the external potential w cx t (r) is such that the clas- 
sical dynamics of a single electron in v ext (r) is integrable. 
The statistics of N non-interacting electrons moving in 
v cx t (r) follows that of the one-particle system, i.e. Poisso- 
nian in this integrable case. Suppose also that the exact 
functionals _E xc [n] and / xc [n](r, r'; to) are known. The 
exact TDDFT procedure must transform the Poissonian 
non-interacting statistics to the correct interacting ones. 
In what follows, when the interacting system displays 
Wigner-Dyson statistics, we refer to this as the P— *WD 
transition. Moreover, in studying how it does so, the 
mechanism for interaction-induced chaos can be better 
understood. 

First, the ground-state KS potential v s (r) is found, and 
we ask: Is the classical dynamics of a single electron in 
v s (r) already chaotic? In this scenario, at least part of the 
interaction-induced chaos appears in DFT as a "chaotic 
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kink" in the Hartree and/or exchange-correlation pieces 
of the potential. The "kink" is that piece of the KS po- 
tential due to which the single-particle classical dynamics 
is chaotic. Excitations from the bare KS potential will 
generally show in this case some degree of level repulsion, 
perhaps enough to agree with the experimental NNS dis- 
tribution. We refer to this as type-KS chaos. 

(We note here that it is preferable to perform the level 
statistics in such a way that N non-interacting electrons 
in a given potential follow the same type of distribution 
as a single electron in that potential: this can be achieved 
if the statistics are performed on energy levels of a fixed 
symmetry class. That is, if there are constants of the mo- 
tion in addition to the total energy, one fixes the value of 
each of those constants. In the case of N non-interacting 
electrons, each orbital eigenvalue is a conserved quantity; 
therefore one fixes all except one, since the sum of the 
orbital energies gives the total energy. This amounts to 
considering only single excitations of the system; these 
can be out of any of the occupied orbitals, depending on 
which were chosen fixed but the resulting statistics will 
be independent of the choice). 

In the second step of TDDFT, the bare KS excitations 
are corrected to the true ones using / HX c : What is its 
effect on the statistics? Certainly if the KS system turned 
out to be integrable, then the entire job of transforming 
the statistics is done by the exact / HXC - We refer to this 
as type-/ HX c chaos. 

For N electrons in one dimension, type-KS chaos is 
impossible. The dynamics of a single electron in one 
dimension is obviously always integrable, since there is 
one conserved quantity (the energy of the electron) and 
one degree of freedom. Using the terminology above, the 
Kohn-Sham potential will never acquire a "chaotic kink" 
in one dimension. In more than one dimension the mod- 
ulation of the external potential that is provided by the 
Hartree and xc terms (Eq. ([3])) may induce chaotic mo- 
tion on the one-electron dynamics. Nevertheless, as we 
will discuss in Section |lVl type-/ HX c chaos is likely to be 
more common also in more than one dimension. Gener- 
ally, most of the job required for the WD transition 
has to be done by / HX c- 

Of course in practise, approximate functionals must 
be used. The question then arises: what properties 
are needed in approximate kernels to achieve the cor- 
rect statistics? Are the present-day approximations good 
enough to capture the P^WD transition as the interac- 
tion is turned on? In particular, the adiabatic approx- 
imation for the kernel: this mixes KS single-excitations 
but cannot fold in multiple excitations [52, IHj]. This 
immediately raises a red flag as a chunk of the true exci- 
tations are excluded from consideration: how significant 
is this chunk? Even if the adiabatic kernel creates level 
repulsion by mixing singles only, if the double-excitations 
compose a significant proportion of the spectra in a cer- 
tain range, then the ATDDFT spectra would not be rep- 
resentative of the true system. The answers to these 
questions are likely system-dependent. We illustrate in 



the next Section most of the concepts discussed so far 
using a model of a one-dimensional quantum dot. 



III. A MODEL 1-D QUANTUM DOT 

Consider the problem of two electrons interacting via a 
soft-Coulomb potential v ee (xi, x%) = (1 + (x± — X2) 2 )~ 1 ^ 2 
in a one-dimensional box of length L. A simple anal- 
ysis shows that the matrix elements of v ee scale ap- 
proximately as 1/L (basis of non-interacting electrons), 
whereas the kinetic-energy matrix elements scale as 1/L 2 . 
The interaction in this sense becomes more important as 
the length of the box is increased. The non-interacting 
system is of course integrable, having two constants of the 
motion (the energies of the individual particles) for two 
degrees of freedom. The weakly-interacting case (small 
L) is then almost-integrable, and the strongly-interacting 
case (large-L) non-integrable. In fact, ref.[25| shows clear 
signatures of interaction-induced chaos in a similar sys- 
tem. We start, however, in the weakly-interacting limit, 
realized for L = 1 a.u. 

The light-shaded bars of Fig.fl]) correspond to the 
nearest-neighbor-spacing (NNS) distribution P(s) ob- 
tained via the following three steps: (1) Get the first 
JV max levels of the energy spectrum {E{\ by exact di- 
agonalization of the 2-electron Hamiltonian, and dis- 
card the lowest -/V min . We typically used N max = 2000, 
^Vmin = 200, and only states of even symmetry; (2) Un- 
fold the staircase function AT(E), representing the num- 
ber of states having energy less than E, by applying the 
map Xi = M{Ei), where Af is the smooth part of J\f, 
found via a 4 th -order polynomial fit to AT(E); (3) Set 
si = (xi + i — x^, and count -P(s), the number of occur- 
rences for each s,/s, where s = J^i s i/( N mnx ~ N min ). 
P(s) then contains information about inherent fluctua- 
tions of the level spacings. As discussed in the introduc- 
tion, these are generally Poissonian for integrable sys - 
tems, and Wigner-Dyson-like for chaotic systems [39j . 
but may also be non-generic for some highly symmet- 
ric systems like the one we are considering here. The 
exact histogram consists of a series of bars separated by 
gaps. The fluctuations around 77(E) (dotted line on the 
inset of Fig.|T])) are evidently not random; they exhibit 
distinctive patterns, each of which gives rise to a differ- 
ent bar in the histogram. In the non-interacting limit, 
the problem of two particles with the same mass in a 
1-dimensional box is equivalent to that of a single parti- 
cle in a 2-dimensional square box, and the patterns just 
mentioned arise from the square symmetry of the prob- 
lem. These remain in the weakly-interacting regime. The 
histogram would be more Poissonian if the two particles 
had different masses, and even more so if the ratio of 
the two masses were an irrational number, as shown in 
Fig.© (and ref.[40i). 

What does TDDFT predict for the histogram within 
the adiabatic approximation?: a single bar centered at 
s = 1 (dark-shaded bar in FigUJ). We now explain why. 
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FIG. 1: Light-shaded bars: NNS histogram for the exact 
spectrum of 2 interacting electrons in a 1-d box of length 
L — 1 a.u. Dark-shaded bar: ATDDFT spectrum within 
exact-exchange. Top arrow indicates that the bar keeps go- 
ing up. Inset: small segment of staircase function J\f(E). The 
dotted line shows the smooth J\f(E). Arrows point to the only 
two single excitations in this region. 




FIG. 2: The NNS distribution for a particle in a 2-d rectan- 
gular billiard of sides a 2 and b 2 depends on the a/b ratio as 
indicated here. The dashed line corresponds to the Poisson 
distribution. (Shown first in ref.[4(J). 



As discussed in Section UH the Kohn-Sharn scheme 
transforms the original problem of 2 interacting elec- 
trons into that of a single electron moving in the po- 
tential v s (x) = v cxt {x) + v nxc (x), where v CKt (x) is the 
external potential (box of unit length, in this case). 
Being a one-dimensional potential, it cannot exem- 
plify type-KS chaos. The Hartree-exchange-correlation 
potential, v HXC (x), was found here within the exact- 
exchange app roximation 41| , as implemented in the octo- 
pus code [44j. v HXC (x) is in this case simply a small bump 



at the bottom of the box (see solid line in lower panel 
of Fig. ([3])), having practically no influence on the high- 
energy spectrum of the non-interacting problem, except 
for a small shift in all the high- lying energies, irrelevant 
when the analysis of the differences between neighboring 
levels is made. The Kohn-Sham equations (EqJ5]) yield 
a set of orbital energies {e,}, but only the subset of oc- 
cupied orbitals is used to get the ground-state energy of 
the interacting system, as: E = J2i occ e i + ^hxcW - 
/ dxn(x)v HXC (x). The excited-state energies E™ DFT are 
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FIG. 3: Ground state density n(x) and Kohn-Sham potential 
v s (x) for 2 interacting electrons in a ld-box, as its length is 
increased from L = 1 to L = 16. The walls of the boxes at 
±L/2 are not shown. 



obtained within linear-response TDDFT as: 



El 



En 



(7) 



where the squares of < DDFT , Sl n = (w™ DFT ) are the 
eigenvalues of the matrix of EqfJJ 

The /hxc kernel we employed is again that of exact- 
exchange, local in time (proportional to 6(t — t')). Its 
Fourier transform with respect to t — t' is therefore 
frequency- independent. As a consequence, from Eq.(f4|) 
and Refs. [Sill], ATDDFT only yields corrections to KS 
single-excitations. Since the frequencies w™ DFT of Eq.® 
are all mixtures of single excitations, the ATDDFT spec- 
trum at this level contains no double excitations at all (or 
multiple excitations for N > 2). This is a serious prob- 
lem even at the low energies that are of interest for most 
quantum-chemical purposes, but in the present setting, 
the missing of double excitations is simply catastrophic. 
Single excitations are not more than a negligible fraction 
of the high-energy set of levels. The arrows of the inset 
in Fig. ([!]), for example, point to the only two single ex- 
citations found in the lowest-energy segment taken into 
account for our statistical analysis. 

We are only one step away from concluding that in 
the energy range entering the statistical analysis, the 
ATDDFT staircase function is identical to the bare KS 
staircase function: 



Af ATDDFT {E)=AT KS (E) 



(8) 



To see this, we just point out that the expectation value 
of Eq. ([5]) becomes vanishingly small in the adiabatic ap- 
proximation for high values of u> q . It can be seen in Table 
U that for unoccupied orbital indices a > 10, adiabatic 
single-excitations barely differ from the KS ones on the 
scale of energy fluctuations. 
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FIG. 4: NNS histograms for 2 interacting electrons in a ld-box 
of length L — 10 7 . P(s) approahes a Wigner-Dyson distri- 
bution (dashed) for very wide boxes, a signature of chaos in 
the underlying classical dynamics. The Poisson distribution 
is indicated by dotted lines. 



Furthermore, in the weakly-interacting limit, Af KS (E) 
coincides with the staircase function for a single electron 
in the presence of v cx t{x), since, as argued before, the 
effect of v axc (x) is negligible. Af KS (E) has no fluctuations 
at all in this limit, being precisely proportional to \/2E, 
and the unfolding process translates this into a histogram 
showing only s = 1 spacings. 

It may be argued that double excitations can simply be 
added to the ATDDFT spectrum by including the sums 
of e a + £b with a and b running over all the unoccupied or- 
bitals. This seems entirely sensible, and works of course 
in the non- interacting limit, but the adiabatic approxi- 
mation yields no corrections to these states, and Eq.© 
still holds true. 

We now discuss what occurs as we move towards the 
strongly- interacting regime. Fig.((l|) shows the NNS his- 
tograms as the length of the box is increased from L = 10 
to L = 10 4 . Strong level repulsion manifests itself with a 
marked decrease of small spacings. The histogram tends 
to a Wigner-Dyson distribution P(s) = %se~^/^ s for 
very wide boxes. This is a quantum signature of the 
underlying interaction-induced chaos, as implied by the 
classical dynamics results obtained by Fendrik et.al.p25l]. 
showing that apart from small regular regions in phase 
space, the dynamics of the system can become strongly 
chaotic due to electron-electron interactions. But it fol- 
lows from Eq.© that ATDDFT yields no level repulsion. 
Practically all levels entering the statistical analysis are 
double excitations. Even if we decided to artificially add 
these by summing unoccupied KS orbital energies, as dis- 
cussed before, only a mixture of uncorrelated levels would 
result, i.e. a Poisson distribution P(s) = e~ s with no 
level repulsion whatsoever (dotted lines in Fig. This 
important signature of interaction-induced chaos is sim- 
ply not captured by ATDDFT. 

As the length of the box increases, localization of the 



TABLE I: Comparison of low-frequency excitations (in 
Hartrees) between bare Kohn-Sham energy differenes e a — U 
and TDDFT frequencies obtained via Eq.Q with the exact- 
exchange approximation to « X c and an adiabatic aproxima- 
tion to /xc- 



i — » 


a 


KS 


TDDFT 


1 -> 


2 


0.0939 


0.2771 


1 -> 


3 


0.3362 


0.5132 


1 -> 


4 


0.6736 


0.8014 


1 -> 


5 


1.1091 


1.1963 


1 -> 


6 


1.6425 


1.7015 


1 -> 


10 


1.7496 


1.7628 


1 -> 


20 


1.9319 


1.9319 



electrons in opposite extremes of the box starts taking 
place, and convergence of the KS equations becomes 
problematic (we used the Octopus code [44[ for the DFT 
calculations). Fig. ([3]) shows the evolution of the KS po- 
tential and ground-state density as L is increased from 
1 to 16. The problems here are no different than those 
encountered when trying to describe the formation of a 
Wigner crystal [45| or other phenomena in the strongly- 
correlated regime. There is already a hint of localization 
when L = 16 (long-dashed line, upper panel of Fig[3]). 
The KS potential adopts the shape of a double well, with 
w HXC (2;) as the barrier in the middle, growing stronger 
with respect to level spacings as the length of the box 
increases. 



IV. DISCUSSION AND OUTLOOK 

Our example model system was an unfortunate one 
for the usual adiabatic approximations of TDDFT. The 
missing double excitations haunt the problem even in the 
integrable limit, making up the significant fraction of ex- 
citations in the true space and being crucial to bring- 
ing about the chaos. There is intensive development to 
go beyond the adiabatic approximation in TDDFT (eg. 
Refs [43|, Hfl H3, El|), and the progress has been success- 
ful in many cases. It will be interesting to adapt the 
kernel of Ref. [43|, |49| to the case in the current paper. 

We can now justify our claim of Section [TT] that type- 
KS chaos is rare. Even in cases where f H xc( r ) develops 
a "chaotic kink" (see second part of Sec|TT|, level repul- 
sion in the bare Kohn-Sham system will only be observed 
when the relevant excitations that enter the statistical 
analysis are of a single-particle nature. In such cases, 
ATDDFT should perform much better than in the model 
of Section IIIII The spectral properties of nonhydrogenic 
atoms in weak external electric or magnetic fields seem to 
be suitable examples. The joint effect of an ionic core of 
inner electrons (describable by a quantum defect for a se- 
ries of single-excitations) and the external field, has been 
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shown to lead to chaos in some parameter regimes [17| . 
The fact that accurate quantum defects can be obtained 
from bare Kohn-Sham potentials, and that the external 
field appears explicitly in the KS equations, suggests that 
ATDDFT has a good chance to succeed in the description 
of this phenomenon. It will be interesting to compare the 
ATDDFT spectrum with that of the scatterer-perturber 
models in the literature, Refs. [l]J H El HI • TDDFT 
treats all electrons quantum mechanically, and provides 
an ab initio method against which to compare the model 
potentials. 

It should also be noted that the level-spacing statis- 
tics physically affect very high-temperature properties 
of the system. Other measures of chaos such as local- 
ization properties of wavefunctions are more relevant at 
usual lower temperatures. Before TDDFT can be used 
for these properties, one must first solve the TDDFT 
"observable problem" i.e. quantities that are directly 
given by the time-evolving density (or in perturbation 
theory, the density response), such as the excitation en- 
ergies, oscillator strengths, or dipole moments, are sim- 
ply extracted from the Kohn-Sham system. However, 
for properties that rely on other aspects of the interact- 
ing wavefunction, such as its localization in phase-space, 



the information is not readily accessible from the Kohn- 
Sham system. Such cases demand the construction of 
appropriate functionals of the time-dependent density for 
the relevant wavefunction-dependent observable, often a 
challenging task (see, for example, Ref. [50(). 

In summary, we have discussed prospects for using 
TDDFT to study how electron-interaction induces chaos 
in an otherwise integrable system. The scalability of 
TDDFT and its success in describing electron correlation 
in an ever- widening range of problems, make it attractive 
for this. Moreover, in principle, TDDFT yields exact 
level statistics for the interacting electronic system. We 
have discussed how the majority of currently-available 
approximate functionals however cannot capture the cor- 
rect statistics in many cases of interest in quantum chaos, 
because of their inability to describe excited states of 
multiple-excitation character. Therefore challenges lie in 
developing approximate and easily-implementable func- 
tionals for this purpose. 
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